DNA Breathing Dynamics in the Presence of a Terahertz Field 
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We consider the influence of a terahertz field on the breathing dynamics of double-stranded 
DNA. We model the spontaneous formation of spatially localized openings of a damped and driven 
DNA chain, and find that linear instabilities lead to dynamic dimerization, while true local strand 
separations require a threshold amplitude mechanism. Based on our results we argue that a specific 
terahertz radiation exposure may significantly affect the natural dynamics of DNA, and thereby 
influence intricate molecular processes involved in gene expression and DNA replication. 



Spectroscopic techniques in the terahertz range are 
currently emerging as new tools for the investigation of 
biological macromolecules [1] . In spite of the experimen- 
tal difficulties (caused by the water's giant Debye dipole 
moment, which leads to a substantial dielectric relax- 
ational loss in the THz range) a notable amount of re- 
search effort is being devoted to the development of so- 
phisticated THz bio-imaging [2]. Nevertheless, very little 
is known about THz-radiation's influence on biological 
systems, and the mechanisms that govern this influence. 
The possibility that low frequency electromagnetic radi- 
ation may affect genetic material, enzymatic reactions, 
etc. was introduced long ago )3 a , and since then has been 
a subject of constant debate. The energy of such radia- 
tion is too low to directly disrupt any chemical bonds or 
cause electronic transitions. Only a resonance- type in- 
teraction might lead to an appreciable, biological effect. 
In biomolecules such interactions are possible through 
the ubiquitous hydrogen bonds that have energies in the 
THz range. Numerous in vivo and in vitro experiments 
have been conducted to clarify low frequency radiation's 
ability to cause biological effects, such as chromosomal 
aberration, genetic damage etc. The experimental stud- 
ies have been conducted under a variety of conditions, 
but mostly at frequencies below 0.01 THz, power below 
1 mW/cm 2 , and short exposure times. The data col- 
lected in these conditions led to mixed conclusions: some 
studies reported significant genetic damages while oth- 
ers, although similar, showed none [4]. The major inter- 
national research project, "THz-bridge" [5], which was 
specifically concerned with THz radiation genotoxicity 
concluded that: under some specific conditions of expo- 
sure, change in membrane permeability of liposomes was 
detected and an induction of genotoxicity was observed 
to occur in lymphocytes. Hence, this project confirmed 
the existence of THz genotoxicity, but it remains unclear 
under which specific conditions such effects occur. 

Recent measurements confirm that only extended (6 
hours) exposure to a weak THz field can cause genomic 
instability in human lymphocytes 6 . Independently, it 
was reported that neurons exposed in vitro to powerful 
THz radiation (over 30 mW/cm 2 ) cause infringement of 



the morphology of the cellular membranes and intracel- 
lular structures [7]. The same work also showed that at 
decreasing power and/or at different frequencies the mor- 
phological changes do not occur. Recently, it was further 
pointed out that exposure to a low level THz radiation 
can interfere with the protein-recognition processes [5]. 
It was also shown that exposure of mice to 3.6 THz high 
power (15 mW/ cm 2 ) radiation, for 30 min caused behav- 
ioral changes [3], while under short (5 min) exposure the 
changes could not be detected. Thus the available ex- 
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FIG. 1: Spontaneous formation of a localized permanent 
opening in dsDNA in the presence of a THz field with ampli- 
tude A = 144 pN, and = 2 and a spatial fluctuation with 
high amplitude (see text). 

perimental data strongly suggest that THz-radiation can 
affect biological function, but only under specific condi- 
tions, viz. high power, or/and extended exposure, or/and 
specific THz frequency. 

The present work is devoted to developing a qualita- 
tive description of the existence of such THz effects. We 
argue that the appropriate conditions correspond to ir- 
radiation parameters that assure the existence of linear 
and nonlinear resonances between DNA conformational 
dynamics and THz field Fig[T] 

We will base our analysis on the Peyrard-Bishop- 
Dauxois (PBD) model of dsDNA [10 , which is arguably 
the most successful available model for describing the 
local DNA pairing/unpairing (breathing) dynamics 
Of particular importance for the THz effect is our previ- 
ous demonstration of strong correlations between regula- 
tory activity, such as protein-DNA binding and transcrip- 
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tion, and the equilibrium propensity of dsDNA for local 
strand separation [12] [T3] HH [15] . With this background 
it is natural to consider the influence of THz radiation on 
dsDNA dynamics within the PBD modeling framework. 
One complication is that the specific physical nature of 
the interactions between DNA and the THz electromag- 
netic field is not known in detail. However, given the 
sensitivity of THz radiation to the strand pairing state 
of DNA [15], we will here simply augment the PBD to 
include a drive in the THz frequency range, without spec- 
ifying the precise nature of the underlying physical cou- 
pling. Since the wavelength of the THz field is larger 
than 100 \xm and the characteristic size the investigated 
dsDNA sequences is less than 20 nm it is reasonable to 
consider an uniform THz field. In this fashion the breath- 
ing dynamics (represented by the normalized separation 
y n of complementary bases) of the n'th base pair in the 
presence of a monochromatic spatially homogeneous ex- 
ternal field is described by 

my n = - U'{y n ) - W'(y n+1 ,y n ) - W'(y n ,y n -i) 

— mjy n + A cos (fit), (1) 

where, the Morse potential U(y n ) = D„(exp (— a n y n ) — 
l) 2 represents the hydrogen bonding of the comple- 
mentary bases. The parameters D n and a n depend 
on the type of the base pair (A-T or G-C). Simi- 
larly, W(y n ,y n -!) = lx[yn,yn-i]{Vn ~ Vn-i) 2 with 
xbjniUn-i] = k(l + pexp(-/3(y„ - y n -i))) represents 
the stacking energy between consecutive base pairs. The 
term mrfy n is the drag caused by the solvent while 
AcosiQt) is the (THz) drive. For simplicity, we con- 
sider here a homogeneous poly(A) [17 DNA molecule 
with 64 base pairs. Our numerical simulations of this set 
of nonlinear coupled equations (using periodic bound- 
ary conditions) showed that the primary response of the 
system is spatially homogeneous (FigJTJthe first 40 ps). 
However, we observed that a presence of fluctuations may 
lead to the creation of persistent spatially localized open- 
ings of the double-stranded molecule, Fig[l] Because, 
such openings (bubbles) are known to functionally affect 
dsDNA we focus here on characterizing the conditions 
under which these states appear. 

Since, the primary response of of the driven system is 
spatially uniform it can be understood in terms of a single 
classical damped and driven Morse oscillator. Such sys- 
tems have been studied quite extensively in many other 
contexts |18j . and are known to display rich dynamical 
behavior including limit cycles, period-doubling bifurca- 
tions etc. An illustration of this is shown in Fig|2j A, 
for our system parameters, in terms of the Poincare sec- 
tion y n (t — mT), where T — 2ir/il, versus the driv- 
ing frequency f2 for several integer values of m. Fcigcn- 
baum's period-doubling route to chaos is clearly seen in 
this figure, which demonstrates that this system's main 
instability mechanism is through period-doubling. As 



unstable 




stable 

S3 tTMzl 



FIG. 2: A) Bifurcation diagram for a poly(A) molecule in 
an external THz field with amplitude A = 144 pN. B) The 
instability curves for wave vector k = (solid curve) and 
k = tv (dashed curve) are plotted. On the vertical axis is the 
amplitude of the external field, and on the horizontal axis is 
the frequency Q of the THz field. 



described earlier, our interest is to investigate whether 
the presences of a THz field can lead to the creation of 
spatially localized unbinding of the DNA double-strand. 
Such states are most easily created in connection with the 
period-doubling events, as a result of a spatially inhomo- 
geneous perturbation. In order to investigate this phe- 
nomenon we adopt a method developed in Ref. [IS]. In 
this work the linear stability analysis (based on the Flo- 
quet theorem) of a single Morse oscillator was performed, 
and it was shown that in a rotating-wave approximation 
|20j the dynamics ensuing from a given perturbation is 
governed by the extended Hill equation. The stability 
analysis of the Hill equation has been developed by Ince 
based on his extension of the Whittaker's method 
of solving Mathieu's equation [35]. Since, our system is 
spatially extended, a perturbation may be spatially in- 
homogeneous. In this case the response of the system 
is characterized by y n (i) = y°(t) + z n (t), where y°(t) is 
the spatially homogeneous and temporally periodic solu- 
tion described above (see Fig. Ill), A, and z n (t) represents 
the perturbations in whose dynamics we are interested. 
In the most general form the perturbation is given in 
terms of Fourier modes: z n (t) — J\ exp[ikn]£k(t) with 
wavenumber, k, and amplitude As with the single 

oscillator, the dynamics of the Fourier amplitudes £fc(i) 
of the perturbation are governed by set of uncoupled ex- 
tended Hill equations. However, the coefficients of these 
equations depend on the wave-number k. The Ince sta- 
bility analysis allows us to analytically determine 23J the 
Floquet exponents (i.e. the stability criteria) as functions 
of the spatial profile as given by k. For more technical de- 
tails on the method see Ref. [24] . The result of our anal- 
ysis are the stability curves dividing the parameter space 
(A/m,n,k) into stability and instability zones. Figure 
|2j panel B shows the stability boundaries of a spatially 
homogeneous perturbation (fc = 0, solid curve) and of a 
staggered (fc = w, dashed curve) perturbation. Perturba- 
tions represented by other wave-numbers have stability 
boundaries in-between those shown in Fig. [2j B. Below 
these curves the spatially homogeneous states perform 
stable temporal oscillations at the frequency of the exter- 
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FIG. 3: The k = tt spatial pattern, at Q = 2.52 THz, is 
presented. The upper panel shows evolution of the system. 
The lower panel shows the power spectrum of the breather 
motion. 
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FIG. 4: Breather formation at Q = 2.00 THz. The upper 
panel shows the evolution of the breather. The lower panel 
shows the power spectrum of the breather motion. 



nal drive and its higher harmonics, while above the curve 
the perturbation leads to period-doubled oscillations. It 
is clearly seen that for driving frequencies below £1 ~ 2.5 
THz spatially homogeneous (k = 0) perturbations lead 
to the strongest instability (occurring at the lowest value 
of the driving amplitude). However, above ~ 2.5 THz 
the staggered (k — tt) perturbation becomes the primary 
instability mode. In the region above the stability bound- 
ary the period-doubled state will therefore tend to adopt 
the spatial profile of the most unstable Fourier-mode (ho- 
mogeneous or staggered). The presence and location of 
the obtained instability curves (Fig. [2] B) were ver- 
ified by exploring the parameter-space [A/m, fi, k) by 
direct simulations of Eq. (1). The dynamics of small 



spatially random (uniformly distributed in the interval 
[— e, e] and thus containing all wave-numbers) perturba- 
tion z n introduced on top of the spatially homogeneous 
and temporally periodic solution y°(t) was followed nu- 
merically. As this perturbation injects energy into all 
wave- numbers, the ensuing dynamics will emphasize the 
mode with the strongest instability. For point A (Fig. 
[2] B), the results are shown in Figure [3j The upper 
panel shows that after the introduction of the pertur- 
bation (at 40 ps) a spatially staggered structure forms. 
The power spectrum of this dynamics (lower panel in Fig. 
[3| , shows that the frequency of the vibrations is mainly 
fi/2, i.e. a period-doubling transition has occurred as 
expected from the above analysis. The power spectrum 
(lower panel) shows that the signature of the drive per- 
sists at f2, in addition to a component at zero frequency, 
which results from the asymmetry of the base-pair po- 
tential U(y n ). A similar scenario is observed for driving 
frequencies below ~ 2.5 THz, but in these cases the 
solutions remain spatially homogeneous after the intro- 
duction of the perturbation. The fact that the instability 
curves for < k < tt lie between the two curves in Fig. [2] 
B shows that it is, in general, impossible for any spatially 
localized feature with a finite (2 base pairs < A < 00) 
length scale A to emerge as a result of a linear insta- 
bility. This is in contrast with many other nonlinear 
condensed matter models, where spatially localized and 
temporally periodic modes (intrinsic localized modes or 
breathers) commonly arise through linear modulational 
instabilities. For our model of dsDNA a perturbation 
of finite amplitude e is required to create localized un- 
bindings of the DNA double strand through nonlinear 
mechanisms. An illustration of this is given in Fig. [4] 
for — 2.0 THz, where the dynamics caused by a ran- 
dom perturbation with e > 0.14 is shown to result in 
a localized opening of the DNA double strand. Again 
the power spectrum, given in the lower panel of Fig. [4] 
shows, as expected, that the localized state vibrates at 
the frequency f2/2 in contrast to the uniform background, 
which remains synchronized to the external drive at f2. 
A more detailed study of the nonlinear mechanisms un- 
derlying the creation of such localized unbinding states 
reveals that a localized injection of a certain amount of 
energy is required. Due to the asymmetry of the base 
pair potential, an effective way to achieve this is a local- 
ized compression of the double strand. Using this insight 
we found the perturbation z n = 0.42 cos (f (n — no)), for 
—4 < n — n < 4 (z n — otherwise) to be an effec- 
tive perturbation for the creation of a localized unbind- 
ing state at no, Figjl] The phase portrait (not shown) of 
the base pair no, at the transient time period just before 
the breather is established, reveals the breather forma- 
tion mechanism. First, the locally injected energy from 
the compression must be sufficient for a few consecutive 
base pairs to undergo a temporary local melting transi- 
tion, i.e. reach the plateau of the Morse potential. Sub- 
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sequently, the stacking interactions and friction dissipate 
the energy to gradually reach a steady state. At given 
driving amplitude and frequency, the dsDNA molecule 
has two steady conformational states with different aver- 
age energies and amplitudes. Only in these states is the 
externally pumped energy exactly compensated by the 
energy dissipated by the friction. The first state is vi- 
brating predominantly at the driving frequency f2, while 
the second mainly oscillates at f2/2. Hence, the breather 
formation proceeds in three stages: 1) all the base pairs 
are in the first state then; 2) a few of them undergo a 
temporary melting transition as a result of the pertur- 
bational energy injection. Finally, 3) the perturbed base 
pairs reach the second steady conformational state, after 
dissipating their additional energy. 

In summary, we have found that a THz field may cause 
dynamical separations of the DNA double-strand. In 
the presence of weak perturbations e.g. thermal fluc- 
tuations, small amplitude response occurs (at half the 
driving frequency) either in a spatially uniform manner 
or, at higher-frequencies, in an unusual spatially dimer- 
ized form. In the latter case, neighboring base pairs os- 
cillate in an out-of-phase fashion. However, large local- 
ized openings (bubbles) in the DNA double strand can 
only occur via a nonlinear mechanism requiring a spatial 
perturbation above a certain amplitude threshold that 
is determined by the intensity and the frequency of the 
THz field. 

We previously showed experimentally that the intro- 
duction of an artificial permanent bubble in dsDNA, via 
a short mismatched segment, is in fact sufficient for tran- 
scription even in the absence of any auxiliary factors 
[15j . Our present finding of large localized openings re- 
sulting from an external THz field therefore underscores 
the importance of including the interactions of genomic 
DNA with the surrounding environment in DNA mod- 
els. The amplitudes of the resulting openings observed 
in the presence of a THz field are significant compared to 
the dynamic signatures of protein binding sites and tran- 
scription start sites [Hj. This suggests that THz radia- 
tion may significantly interfere with the naturally occur- 
ring local strand separation dynamics of double-stranded 
DNA, and consequently, with DNA function. 

Based on the model results present here, we believe 
that the main effect of THz radiation is to resonantly 
influence the dynamical stability of the dsDNA system. 
Hence, our instability curves define the parameters space 
(i.e. the amplitude of the field (or the power), and fre- 
quency) in which THz radiation can have an immediate 
effect. In contrast, nonlinear instability may occur at any 
point of the instability diagram but it requires fluctua- 
tions with significant amplitudes. In biological systems, 
such fluctuations arc generated thermally. Hence, the oc- 
currence of a fluctuation with sufficiently large amplitude 
is very rare, and therefore extended exposure is required 



for the THz effect to take place, via nonlinear instabili- 
ties, especially if the power is small. In this framework, 
it is than natural that the character of THz genotoxic 
effects are probabilistic rather than deterministic. 
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